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The collapse of slowly rotating molecular cloud cores threaded by magnetic 
fields is investigated by high-resolution numerical simulation. Outflow formation 
Oh 1 in the collapsing cloud cores is also followed. In the models examined, the cloud 

core and parent cloud rotate rigidly and are initially threaded by a uniform mag- 
netic field. The simulations show that the cloud core collapses along the magnetic 
field lines. The magnetic field in the dense region of the cloud core rotates faster 
than that of the parent cloud as a consequence of spin-up of the central region 
^ ■ during the collapse. The cloud core exhibits significant precession of the rotation 

axis, magnetic field, and disk orientation, with precession highest in the models 
with low initial field strength (< 20 /xG). Precession in models with initial fields 
of ~ 40 (jlG is suppressed by strong magnetic braking. Magnetic braking trans- 
fers angular momentum form the central region and acts more strongly on the 
component of angular momentum oriented perpendicular to the magnetic field. 
After the formation of an adiabatic core, outflow is ejected along the local mag- 
netic field lines. Strong magnetic braking associated with the outflow causes the 
direction of angular momentum to converge with that of the local magnetic field, 
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resulting in the convergence of the local magnetic field, angular momentum, out- 
flow, and disk orientation by the outflow formation phase. The magnetic field of 
a young star is inclined at an angle of no more than 30° from that of the parent 
cloud at initial field strengths of ~ 20 fiG, while at an initial field strength of 
~ 40 fiG, the magnetic field of the young star is well aligned with that of the 
parent cloud. 

Subject headings: MHD — methods: numerical — stars: formation — ISM: jets 
and outflows — ISM: clouds 



1. Introduction 

Magnetic fields are believed to play an important role in the gravitational collapse of 
molecular cloud cores. The characteristics of the outflow and jet associated with young 
stars are related to the magnetic field and the rotation of the (proto)star and circumstellar 
disk. Observations of polarization have suggested that the outflows and jets often tend to 
be aligned with each other and with the cloud-scale magnetic field (e.g., Cohen, Rowland, & 
Blair 1984; Strom et al. 1986; Vrba et al. 1986; Vrba, Strom, & Strom 1988; Tamura & Sato 
1989; Jones & Amini 2003). The polarization of young stellar objects (YSOs) also suggest 
that circumstellar dust disks around young stars are aligned perpendicular to the magnetic 
field (e.g., Moneti et al. 1984; Tamura & Sato 1989). Recent high-resolution observations 
of submillimeter polarization have resolved the magnetic fields around young stars on the 
~ 1000 AU scale, which is comparable to the outflow scale, and interactions between the 
magnetic fields, the outflows, and the envelopes have been inferred (Momose et al. 2001; 
Henning et al. 2001; Wolf, Launhardt, & Henning 2003; Vallee, Greaves, & Fiege 2003). 
Wolf et al. (2003) suggested that two Bok globules are associated with outflows parallel to 
the magnetic field, while the two other globules are associated with outflows perpendicular 
to the magnetic field. Thus, the findings for the correlation between outflow, envelope, and 
magnetic field have yet to converge. 

Although the alignment of outflows and jets with the magnetic field has been discussed 
intensively in observational studies, very little theoretical consideration has been given to 
this problem. Theoretical studies of the collapse of a magnetized cloud have been restricted 
to models of aligned rotators, in which the rotation axis is parallel to the magnetic field (e.g., 
Nakamura, Hanawa, & Nakano 1995; Tomisaka 1996; Nakamura, Matsumoto, Hanawa, & 
Tomisaka 1999). Several other studies (Mouschovias & Paleologou 1979; Dorfi 1982, 1989) 
have investigated inclined rotators. However, a collapsing inclined rotator has not been 
investigated for at least 15 years. 
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It is believed that ambipolar diffusion may play an important role in the gravitational 
collapse of cloud core. Theoretical studies including the ambipolar diffusion have done by 
Mouschovias and the co-authors (e.g., Fiedler & Mouschovias 1993; Ciolek & Mouschovias 
1994), and recently by Nakamura & Li (2003). However, they simulated collapse of aligned 
rotators by two-dimensional simulations. The three-dimensional simulations are necessary 
to investigate alignment of an outflow with magnetic field in a collapsing inclined rotator. 
At this moment, three-dimensional simulations including ambipolar diffusion just begin by 
Hosking & Whitworth (2004) using a SPH scheme. In order to investigate collapse of inclined 
rotators, we perform three-dimensional simulations assuming the ideal MHD. Nakano, Nishi, 
& Umebayashi (2002) shows that timescale of magnetic flux loss is considerably longer than 
the free-fall time in n < 10 11 cm -3 , indicating that the ideal MHD is valid during the early 
phase of dynamical collapse. 

Observations suggest that the rotation axis of the molecular cloud is not generally 
aligned with the magnetic field, as exemplified by the filamentary molecular cloud of Taurus 
molecular cloud 1 (TMC1), where the magnetic field is perpendicular to the filamentary axis 
(Tamura & Sato 1989) and a velocity gradient perpendicular to the filamentary axis has 
also been observed for TMC1 (Olano, Walmsley, & Wilson 1988). Accordingly, molecular 
cloud cores, as substructures of molecular clouds, probably have inclined rotational axes 
with respect to the cloud-scale magnetic field. Goodman et al. (1993) indicates that the 
orientation of rotation axes of cloud cores distribute randomly in the Taurus region, and 
suggests that the rotation axes are uncorrelated with the cloud-scale magnetic field. 

Simulations of inclined rotators have also been performed, but with the objective of 
studying magnetic braking of molecular clouds rather than molecular cloud cores (Dorfi 
1982, 1989). In these simulations, the magnetic field is artificially anchored to the boundary 
of the computation box, and only the cloud of interest rotates, while the surrounding gas 
rests at the initial state. The latter assumption induces strong shear flow at the cloud 
boundary, which promotes magnetic braking. Furthermore, the initial state also requires 
some explanation. If a cloud core forms from a turbulent cloud, both the cloud core and 
the envelope must rotate. In the present study, realistic models in which both the cloud 
core and parent cloud rotate and which are initially threaded by a uniform magnetic field 
are considered. The collapse of the cloud core in these models produces differential rotation 
even from a uniformly rotating core-envelope system, which promotes magnetic braking. 

Theoretical studies of outflow formation in the collapsing core have also been restricted 
to axisymmetric models (Tomisaka 1998, 2002; Allen, Li, & Shu 2003). Axisymmetric mod- 
els are applicable only if the rotation axis and magnetic field converge to be parallel, related 
to magnetic braking before outflow formation. As is shown by the simulations in this study, 
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convergence cannot be completed before outflow formation. The effects of outflow on con- 
vergence are considered here in some detail. 

In this paper, the gravitational collapse of slowly rotating magnetized cloud cores 
threaded by inclined magnetic fields is investigated. In § 2 and § 3, the basic equations 
governing the simulations and the model of magnetized molecular cloud cores are intro- 
duced. In § 4, the methods of numerical simulations are presented, and in § 5, the results 
of the simulations are shown and the inclination of the magnetic field, angular momentum, 
disk, and outflows are estimated. In § 6, precession and alignment of several axes of a YSO 
are discussed with respect to magnetic braking, and outflow velocities and momenta are 
considered. Finally, the paper is concluded in § 7. 



2. Basic Equations 

Self-gravitational magnetohydrodynamic (MHD) equations are assumed, as described 

by 

^ + V-(pv) = 0, (1) 
^( P «) + V-( P ™ T ) + V(P + ^ -V- (^-) +PVM> = 0, (2) 
^B + V • (vB T - Bv T ) = , (3) 

and 

V 2 ^ = AnGp , (4) 

where p, v, P, B, and G denote density, velocity, pressure, magnetic field, gravitational 
potential, and the gravitational constant. The transpose of B is denoted by B T . The 
equation of state for a molecular cloud is well described by an isothermal gas for low-density 
regions and a polytropic gas for high-density regions (Tohline 1982). Two types of equation 
of state are therefore employed, i.e., 

P= { ° 2 sP for P < / 5 x 

I c 2 sPcr (p/p cr ) 5 / 3 for p > p cr ' 

where c s and p cr denote the isothermal sound speed and the critical density. The gas tempera- 
ture is assumed to be 10 K (c s = 0.19 km s" 1 ) below the critical density p cr = 2 x 10~ 13 g cm -3 
(n cr = 5.24 x 10 10 cm -3 ), and to increase adiabatically above it. 
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3. Models of Cloud Cores 

As an initial model of a molecular cloud core, we consider a slowly rotating, spherical, 
isothermal cloud threaded by a uniform magnetic field. The cloud is confined by a uniform 
ambient gas that mimics a parent molecular cloud. Both the cloud and the ambient gas have 
a uniform magnetic field and rotate with the same angular rotation speed. 

As a template for a molecular cloud core, we consider the density profile of the critical 
Bonnor-Ebert sphere (Ebert 1955; Bonnor 1956). If £>be(£) denotes the non-dimensional 
density profile of the critical Bonnor-Ebert sphere (see, Chandrasekhar 1939), the initial 
density distribution is given by 

( ) = I Po0BE(r/a) for r < R c 
P[T) ~ I PoQBE(R c /a) for r>R c ' W 



and 



, x 1/2 

' , (7) 



where r, /, and po denote the radius, density enhancement factor, the initial central density, 
respectively. The critical Bonnor-Ebert sphere is obtained when f — 1. The radius of the 
cloud is defined by R c = 6.45a, where the numerical factor comes from the non-dimensional 
radius of the critical Bonnor-Ebert sphere. The density contrast of the initial cloud is 
p(0)/p(R c ) = 14. 

In this paper, p is set equal to 1 x 1CT 19 g cm~ 3 , which corresponds to a number density 
of n = 2.61 x 10 4 cm~ 3 for the assumed mean molecular weight of 2.3. The initial freefall 
timescale at the center of the cloud is thus tg = (3n /32Gp ) 1 ^ 2 = 2.10 x 10 5 yr. The density 
enhancement factor is set at / = 1.68 to allow the rotating magnetized cloud to collapse. 
The radius of the cloud is thus R c = 0.178 pc, and the mass of the cloud is 6.130 M & . The 
parameter a (= -Eth/I-Egravl) is equal to 0.5, where E th is the thermal energy of the cloud 
and .Egrav is the gravitational energy. 

The initial cloud and ambient gas rotate rigidly with angular velocity of Qq and rotation 
axis parallel to the z axis. We assume Qo — 7.11 x 10 _7 yr _1 (= l.SOt^ 1 ), which corresponds 
to f3 = 0.02, where (5q = £ , rot /|£ l grav | and E TOt is the rotation energy. This value of f3 is 
typical for rotating molecular cloud cores (Goodman et al. 1993). 

The initial magnetic field is given by 

/ sin6> \ 

(8) 



B x \ 


( sin# 
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\= B \ 


B z J 


\ cos 6* 
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in Cartesian coordinates (x,y,z), where B denotes the initial magnetic field strength and 
9 denotes the angle between the initial magnetic field and the initial angular momentum. 
When 9 = 0, the magnetic field is aligned with the angular momentum. Various models can 
be constructed by changing the model parameters B and 9, allowing the effect of strength 
and configuration of the magnetic field on the collapsing cloud cores to be investigated. 

Models with four different magnetic field strengths were constructed, as shown in Ta- 
ble 1. The model names are derived from the initial strength of the magnetic field, strong 
(SF), moderate (MF) or weak (WF), and the initial angle 9 (00-90), with the no- field 
model denoted NM. The initial strength of a magnetic field is characterized by three non- 
dimensional quantities, B /B cr , -Emag/I-Egravl, and (3 P . The first quantity B /B cr is the ratio 
of field strength to critical field strength, where the critical field strength is defined by 
B CT = G 1 / 2 S/27T and £ is the central column density of the cloud given by £ = 2 p(r)dr. 
When the initial magnetic field is weaker than the critical magnetic field, the cloud is mag- 
netically super-critical, leading to collapse as the gravity exceeds the magnetic pressure 
and tension without sufficient thermal pressure support. This paper considers the case of 
super-critical clouds. The second quantity -E m ag/|-Egrav| is the ratio of magnetic energy to 
gravitational energy, where E mag = J B^/ft-ndr. Equipartition of magnetic energy and grav- 
itational energy is achieved (E mag /\E gTav \ = 1) when B = 21.9 /xG in the present model. 
All models except for the weak field models (B = 7.42 /xG) and model NM have magnetic 
energy comparable to the gravitational energy. The third quantity (3 P is the plasma beta at 
the center of the cloud, as defined by f3 p = 8ttc 2 s p /Bq, which takes values of the order of 
unity except in the weak field models and model NM. 



4. Numerical Methods 

A self-gravitating MHD simulation code was developed by extending the simulation 
code of Matsumoto & Hanawa (2003b). The MHD equations and Poisson equation are 
solved by a finite difference method with second-order accuracy in space and time. A nested 
grid is employed to solve the central region with higher spatial resolution. The nested grid 
consists of a concentric hierarchical Cartesian grid (Matsumoto & Hanawa 2003b; Machida, 
Tomisaka, & Matsumoto 2004), and the cell width of each grid decreases successively by a 
factor of 2. The coarsest grid is labeled level / = 1, and the Ith level grid has 2' -1 times 
higher spatial resolution than the coarsest grid. The sub-grids are automatically formed 
when Aj/16 becomes smaller than the cell width of the finest grid, where Aj denotes the 
Jeans length (c.f., Jeans condition of Truelove et al. 1997). The maximum level of the grid is 
^max = 14 in typical models. The detail implementation is shown in Matsumoto et al. (2004). 
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The numerical fluxes of the MHD equations (1) - (3) are obtained using the linearized 
Riemann solver reported in Fukuda & Hanawa (2000). This solver is based on the method in 
Roe (1981) but has been extended for solving ideal, barotropic MHD equations. A MUSCL 
approach and predictor-corrector method are adopted here for integration over time in order 
to achieve second-order accuracy in space and time (e.g., Hirsch 1990). Similar to standard 
adaptive mesh refinement (AMR), a multi-timestep scheme is adopted here (e.g., Chiang, 
van Leer, & Powell 1992). When the calculation of a fine grid is synchronized with that for 
coarser grids, the numerical fluxes of the coarse grid at the grid interfaces (between fine and 
coarse grids) are modified by the fluxes calculated by the fine grid in order to conserve the 
numerical fluxes at the grid interfaces. In order to keep V B = within a round-off error, the 
divergence-free staggered mesh of Balsara (2001) and Balsara & Spicer (1999) is employed. 
This method is similar to the constrained transport approach of Evans & Hawley (1988), 
but has been optimized for the Godunov-type Riemann solver and hierarchical grids. This 
method makes it possible to stably solve even extremely low-(3 regions of f3 p = 10~ 5 — 10~ 3 . 

Self-gravity is updated in every step across all grid levels by solving the Poisson equation 
(4) using the multi-grid iteration for a nested grid developed by Matsumoto & Hanawa 
(2003a). In this code, a full multi-grid scheme is employed (see the algorithm shown in Press 
& Teukolsky 1991) to accelerate the convergence of the red-black Gauss-Seidel iteration. 

To test the MHD code, we computed MHD shock tube problems introduced by Brio 
& Wu (1988), confirming that every wave propergates across the grid interfaces without 
reflection. Several tests for self-gravity are described in Matsumoto & Hanawa (2003a). 

A boundary condition is imposed at r = 8R C , where the magnetic field and ambient gas 
rotate at an angular velocity of Qq. The velocity and magnetic field at the boundary are 
therefore given by 

v b (t) = Vl z x R , (9) 

and 

B b (t) = B cos(lV) + z(z ■ B ) [1 - cos(fi t)] + B x z sin(fi t) , (10) 

where z denotes the unit vector in the z direction and R denotes a position vector of the 
cylindrical radius. The density is fixed whenever equal to the initial density of p(R c ) at this 
radius. The large radius of this boundary condition eliminates the effect of artificial magnetic 
braking due to the boundary condition. Even for a strong field model of Bq = 37.1 fiG, the 
crossing time for an Alfven wave from the center to the boundary is 1.1 x 10 6 yr, which is 
5.3 times longer than the freefall time. Gas is considered to be the only source of self-gravity 
in r < R c when solving the Poisson equation (4). Point-symmetry with respect to the origin 
is imposed in order to reduce computational cost. 
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Each grid consists of 128 x 64 x 128 cubic cells in an (x, y, z) coordinate frame. At the 
initial stage, five sub-grids (/ max = 5) are set, where two fine sub-grids contain the cloud 
and the remaining three (sub-)grids contain the ambient gas. The cloud with radius R c is 
therefore solved with a grid of 4 < / < Z max , and the cloud is resolved by 128 x 64 x 128 x 

(Imax - 3) Cells. 



5. Results 

5.1. Clouds Threaded by Moderately Strong Magnetic Fields 

5.1.1. Model MF 45 

Figure 1 shows the evolution of model MF45 (9 = 45°) as a series of cross sections. As 
shown in Figure 1/, the plane of each cross section is defined so as to include the z and x 
axes for planes A and B, oriented parallel to the mean magnetic field in the central dense 
region of p > 0.1p max , where p max denotes the maximum density. The mean magnetic field 
in the central dense region (the central magnetic field) is defined by 

B c =\\ B(r)dr, (11) 

max 

where V = J p>01p dr denotes the volume of p > 0.1p max . Plane A forms an angle 9 Z 
with the x — z plane, and plane B forms an angle 9 X with the x — y plane, where 9 Z = 
arctan(fi Cj j / /_B Cia; ) and 9 X = axctan(B CtZ / B c ^ y ). 

Figure la shows the structure of the entire cloud core in terms of the cross section 
of plane A at p max = 5.58 x 1CT 18 g cm" 3 . The cloud undergoes runaway collapse, and 
the central density increases. The cloud has an oblate density structure with major axis 
extending almost perpendicular to the local magnetic field, as the gas falls along the field 
lines. The cloud becomes thinner in the denser region, and collapses with supersonic flow in 
the direction perpendicular to the major axis. The central magnetic field of the dense region 
rotates around the z axis by up to 9 Z = 32.7°, while the initial magnetic field is parallel to 
the x — z plane [9 Z = 0). 

Figure lb shows the central cloud at p max = 1.70 x 1CT 14 g cm" 3 . The cloud forms 
a collapsing isothermal disk at the center, oriented perpendicular to the local magnetic 
field. The central magnetic field rotates by up to 9 Z = 37.7° by this stage. This collapsing 
magnetized disk resembles that in axisymmetric simulations by Nakamura et al. (1995) and 
Tomisaka (1996). 
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When the central density exceeds the critical density p cr , an adiabatic core forms at the 
center of the disk. The adiabatic core terminates the runaway collapse, and the accretion 
phase begins, in which the adiabatic core accretes gas from the disk-like infalling envelope, 
as shown in Figure lc. The adiabatic core has a radius of 7 AU and corresponds to the first 
core of Larson (1969). 

The core spins faster than the infalling envelope, causing the magnetic field to become 
twisted. A bipolar outflow then forms approximately 400 yr after formation of the adiabatic 
core. The evolution of this outflow was followed for ~ 600 yr. By the final stage, the 
Alfven velocity increases up to 189 km s _1 (= 945c s ), which restricted the timestep of the 
simulation due the Courant condition and rendered it necessary to halt the simulations at 
this stage. Figures lc? and le show different cross sections at the final stage of the simulation. 
The outflow can be seen to have extended out to ~ 200 AU in the final stage, remaining 
anchored to a small disk of radius r ~ 20 AU around the adiabatic core. The outflow is 
perpendicular to the infalling disk and has an hourglass structure. Gas is ejected from the 
small disk, but gas also continues to fall into the adiabatic core through narrow channels 
near the polar axes. The field lines are considerably twisted in the outflow region, as seen 
in Figures 2c and 2d. The outflow has a slightly asymmetric structure with respect to the 
disk axis. The high-velocity gas has a spiral structure that traces the magnetic field lines, 
resembling that of Tomisaka (1998, 2002) except for the non-axisymmetric structure. The 
outflow in this case is categorized as a U-shaped outflow according to Tomisaka (2002). 

Figure 2 shows three-dimensional views of the cloud core at the final stage at three 
different scales. In Figure 2 a, which shows the entire region of the cloud core at the same 
scale as Figure la, the disk of the cloud core is clearly perpendicular to the magnetic field, 
which is anchored to the outer boundary and is rotated around the z axis by 19°. The 
central disk axis, however, is rotated by ~ 60° due to spin-up of the central region, as seen 
in Figure 2 c, but the disk remains perpendicular to the magnetic field. 

The spin-up of the central cloud can be investigated by evaluating the angular velocity 
of the cloud center as follows. 



\ u c\ 

M(a\ + of) 



(12) 



where the central angular momentum is defined as 




'max 



r x vpdr , 



(13) 



and the mass is defined as 




'max 



pdr . 



(14) 
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Here, a n (n — 1, 2, 3) denote the axis lengths of the dense region, which can be evaluated by 
measuring the specific moment of inertia, i.e., 




hj — ~TJ I r^pdr , (15) 

m </p>0.1p max 

following the analysis of Matsumoto & Hanawa (2003b). The subscripts % and j represent 
coordinate labels, i.e., x — n, y — r 2 , and z = r 3 . The three axis lengths a n are defined as 
the square roots of three eigenvalues A n of the specific moment of inertia 1^, as given by 

a n = XH 2 , (16) 
= P l IP , (17) 

V A 3 J 

and 

P = (Pi,P2,P 3 ) , (18) 

where p n denote the eigenvectors of I. The axis lengths a n correspond to the mass-weighted 
radii of interest. Three eigenvectors p n correspond to the directions of the principal axes 
of the dense region. Arranging the eigenvalues in descending order (Ai > A 2 > A3), the 
minimum eigenvalue a 3 denotes the minimum axis and corresponds to the thickness of the 
disk-like cloud. The associated eigenvector p 3 denotes the normal vector of the disk-like 
cloud. The eigenvector can used for evaluating the inclination of the disk, as shown later. 

Figure 3 shows the evolution of the central angular velocity f2 c defined above and the 
central magnetic field strength \B C \ (see eq. [11]). It can be clearly seen that the angular 
velocity and magnetic field strength increase in proportion to pU&x up to the stage of adiabatic 
core formation (see lines for model MF45), implying that the cloud undergoes self-similar 
collapse, maintaining its disk shape (see Nakamura et al. 1995, 1999). In the early accretion 
phase (p max > Per), the growth of the magnetic field and angular velocity are faster than in 
the runaway collapse phase. After outflow formation (p max — 10~ n g cm -1 ), this growth 
slows. 

As collapse proceeds, the central dense region spins up faster than the less-dense region. 
The magnetic field threading the dense region is thus wound up and twisted, causing the 
angular momentum to be redistributed, and it is this process that is called magnetic braking. 
In order to estimate the magnitude of magnetic braking, the angular momentum in the 
dense region was followed more closely. Figure 4 shows the evolution of the central angular 
momentum after decomposition into two components, parallel and perpendicular to the 
magnetic field, as defined by 

J,, = \J C -B C \/\B C \ , (19) 
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and 

Jl = | J c x B C \/\B C \ , (20) 

respectively. 

As collapse proceeds, J±/M 2 decreases considerably, while J\\/M 2 decreases only slightly 
in the runaway collapse phase. The ordinate, J/M 2 , approximated to the specific angular 
momentum per mass, j/M, remains constant as the disk-like cloud collapses in the cylin- 
drically radial direction (see, e.g., Matsumoto, Hanawa, & Nakamura 1997; Matsumoto & 
Hanawa 1999). The decrease in both components is therefore attributed to magnetic braking, 
by which angular momentum is transferred from the central dense region to the surround- 
ing infalling envelope. In the range 10 _14 gcm -3 < p max < 10~ 12 gcm~ 3 , J±/M 2 increases 
due to transfer of the angular momentum back from the infalling envelope to the dense re- 
gion. Similar undulations of angular momentum also occurred in the perpendicular rotator 
of Mouschovias & Paleologou (1979). In the accretion phase, J±/M 2 decreases drastically 
after outflow formation. In the outflow formation phase, the angular momentum is there- 
fore aligned with the magnetic field. The parallel component J\\/M 2 continues to decrease 
considerably in this phase, similar to the behavior in axisymmetric models (e.g., Tomisaka 
1998, 2002). 

Figure 5 shows the loci of the central magnetic field B c , the central angular momentum 
J c , and the normal vector of the disk-like central cloud p 3 . Each vector is plotted in the 
two-dimensional plane as arctan(V 2 /\4 y )\4~ 1 (V r a; , V y ), where V = (V x ,V y ,V z ) represents a 
given vector and V xy = (V 2 + V 2 ) 1 ^ 2 . The distance from the origin represents the angle 
between a given vector and the z axis. For example, a vector parallel to the x axis is plotted 
at (90°, 0°). The direction of the magnetic field at the boundary r = 8R C is also plotted 
to represent the direction of the global magnetic field B\>. This magnetic field rotates at a 
constant angular velocity Q . 

The central magnetic field starts from the point (45°, 0°), and the central angular mo- 
mentum starts from the origin according to the initial condition. The normal direction of 
the density structure is located at ~ (30°, 0°) in the very beginning of the runaway collapse 
phase. At adiabatic core formation, the central magnetic field rotates around the z axis 
with precession of up to 28° in the runaway collapse phase, while the direction of the global 
magnetic field rotates around the z axis by only 19° in this phase. This difference in the 
direction of the central and global magnetic fields is due to the spin-up of the dense region. 
The angular momentum also undergoes precession, gradually coming to be aligned with the 
central magnetic field as J± decreases faster than Jii (Fig. 4). The normal vector of the disk 
is also aligned with the central magnetic field, suggesting that the cloud collapses along the 
magnetic field lines, and the disk-like structure is formed perpendicular to the local magnetic 
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field. At the stage of adiabatic core formation (p max = Per), the angle between the central 
angular momentum and the central magnetic field increases, due to the increase of J± up 
to the stage of p max ~ 1CT 12 cm g -1 (Fig. 4). During outflow formation, B c , J c , and p 3 
begin to come into alignment, and the three vectors eventually become parallel. The central 
magnetic field B c rotates around the z axis by 9 Z = 62°. 

Figure 6 shows the evolution of B c , J c , and p 3 with respect to the global magnetic 
field i?b- In this figure, the origin represents the direction parallel to the global magnetic 
field. In the early stage of the runaway collapse phase, J c and p 3 pursue B c , and thereby 
move toward the origin in this plot. However, B c moves to the right in this plot, reflecting 
the spin-up of the central cloud. In the final stage, the three vectors are inclined from the 
global magnetic field by ~ 30°. As shown later, this angle decreases with increasing initial 
magnetic field strength. 

Figure 7 shows the overall structure of the cloud core in the final stage. Figure 7 a shows 
the directions of the mean magnetic field, mean angular momentum, and normal vector of 
the disk as a function of radius. These vectors are evaluated in the same manner as for 
B c , J c and p 3 , except that the volume of integration is bounded within a given radius in 
equations (11), (13), (14), and (15). The magnetic field and normal vector of the disk are 
almost parallel at all radii, whereas the angular momentum gradually becomes aligned with 
the other directions as r decreases from ~ 10 4 AU to r ~ 100 AU. In r < 100 AU, all three 
vectors are roughly aligned, and the magnetic field and normal vector coincide remarkably. 

It is noteworthy that the curves in Figure 7a are quite similar to those in Figure 5 
except for the lack of precession in the former. This indicates that the magnetic field, 
density distribution, and angular momentum distribution are affected in part by the history 
of cloud collapse. The early distribution is maintained in the envelope during evolution in 
the central dense region, as the less-dense envelope has a longer dynamic (freefall) timescale. 

Figure 7b shows the mean magnetic field strength and rotation velocity with respect 
to radius r. The magnetic field strength increases monotonically up to ~ 1 G toward 
the center, with amplification of up to ~ 10 5 times. In contrast, the rotation velocity 
gradually decreases from the boundary to 10 3 AU, but increases rapidly in the interval 
between 10 3 AU and 10 AU, reaching a maximum of ~ 0.2 km s -1 , which is approximately 
equal to the sound speed c s . The rotation curve for r < 10 AU reflects the rigid rotation of 
the adiabatic core. This rotatin velocity is considerably lower than the Keplerian velocity of 
~ 2 km s _1 at ~ 10 AU, indicating that the adiabatic core and the surrounding disk are not 
rotationally supported. However, the gas in the outflow region spins at the rotation velocity 
of 2 — 4 km s _1 , roughly equal to the outward radial velocity of the outflow. This high-speed 
rotation disappears from Figure 7b in the process of averaging. 



-13- 



5.1.2. Models with Different 9 

The effect of the angle 9 is investigated by comparing models MFOO, MF45, and MF90. 
Model MF90 is an extreme model, in which the magnetic field is aligned perpendicular to 
the angular momentum in the initial stage. The results for MF70 and MF80 are also shown 
for comparison with this extreme model. 

Figure 3 shows the central magnetic field B c and the central angular velocity Q c as a 
function of maximum density p max for models MFOO, MF90 and MF45. For all models, the 
amplification of the magnetic field due to cloud collapse is almost identical. The magnetic 
field increases in proportion to plll x in the runaway collapse phase, irrespective of 9. This 
suggests that the central disk threaded perpendicularly by the magnetic field collapses along 
the major axes (cylindrically radial direction), regardless of 9. However, the evolution of 
the central angular velocity Q c does depend on 9, and models with larger 9 have suppressed 
amplification of Q c . In particular, Q c undulates and changes sign around the stage of p max ~ 
10~ 14 g cm -3 in model MF90 due to strong magnetic braking. In this extreme model, the 
angular momentum remains parallel to the z axis due to the symmetry assumed in the 
initial stage. This change in the sign of Q c represents anti-rotation of the central cloud (c.f., 
Mouschovias & Paleologou 1979). 

Figure 4 shows the evolution of the central angular momentum J c by decomposition 
into parallel and perpendicular components with respect to the central magnetic field. The 
parallel component J\\/M 2 of model MFOO slightly decreases in time, similar to model 
MF45. Model MFOO has only a parallel component, as the perpendicular component van- 
ishes (Jj_ = 0) due to the axisymmetry of the model. The perpendicular component J±/M 2 
of MF90 decreases considerably, also similar to model MF45. After outflow formation, Jj_ is 
expelled almost completely from the central region and adiabatic core. Therefore, the mag- 
netic braking against each component is roughly independent of 9, with the perpendicular 
component being more strongly affected by magnetic braking than the parallel component. 

The total angular momentum of model NM, which has no magnetic field, is also plotted 
for comparison in Figure 4. The angular momentum increases slightly, roughly oc pUax, as 
collapse proceeds. The slight increase in J/M 2 for model NM seems to originate from the 
incorporation of a spherical collapse with the disk-like collapse, resulting in an amplification 
of J/M 2 (see, e.g., Matsumoto et al. 1997; Matsumoto & Hanawa 1999). The slight decrease 
in J\\/M 2 for models MFOO and MF45 is due to weak magnetic braking. 

Figure 8 shows the outflows seen in model MF90. Surprisingly, this extreme model 
also forms an outflow, where the outflow is aligned parallel to the local magnetic field and 
perpendicular to the disk. The outflow has a flat shape that is restricted to the region near 
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the x — y plane, as seen in Figure 8a, while the angular momentum is parallel to the z 
direction. This outflow is quite different from that in models with 9 ^ 90°. In this model, 
the spin of the central adiabatic core trails the threading magnetic field. The field lines thus 
pinch into the region near the equatorial plane and bend in the plane. Consequently, the 
gas is accelerated in the plane. In contrast, outflow is associated with twisted magnetic field 
lines in models with 9 ^ 90° (Figs. 2c and 2d). 

Models MF70 and MF80 were constructed in order to evaluate off-symmetry config- 
urations such as MF90 and MF00. Figure 9 shows the evolution of B c , J c , and p 3 for 
models MF70 and MF80. These progressions are qualitatively similar to that in model 
MF45 (Fig. 5). The normal vector p 3 is more closely aligned with the central magnetic field 
B c in these models with larger 9, and these vectors still drift with precession. The angular 
momentum J c also exhibits precession, but the amplitude of drift is much greater than for 
other two vectors. This precession of J c corresponds to the undulation of J± during collapse. 
After outflow formation, J c inclines to become parallel with the other vectors B c and p 3 . 
In the final stage, the three vectors are all roughly parallel. The convergence points form an 
angle of ~ 35° with the global magnetic fields in both models. 



5.2. Clouds Threaded by Strong Magnetic Fields 

The models with strong magnetic field were initialized with a field strength of 37.1 /iG. 
The evolution of the cloud cores in these models was followed until immediately prior to 
outflow formation. Extremely high Alfven velocities (e.g., 430c s = 82 km s" 1 for model SF00) 
were achieved in the outflow formation stage, which severely restricted the timestep of the 
simulations due to the Courant condition and rendered it necessary to halt the simulations 
before outflow formation. 

Figure 10 shows the central cloud in the final stage for model SF45. The density and 
velocity distributions resemble those of model MF45 (Fig. lc). The density distribution 
reveals a disk shape oriented perpendicular to the local magnetic field, although this is not 
shown explicitly in the figure. 

Figure 11 shows the loci of B c , J c , and p 3 for model SF45. In contrast to model MF45, 
precession is suppressed, and the three vectors are fully aligned by the end of the runaway 
collapse phase, indicating that J± is almost entirely expelled from the region of interest, 
leaving only Jy. Furthermore, the central magnetic field forms an angle of only 5.9° with 
the global magnetic field in the final stage, indicating that the internal magnetic field of the 
cloud core changes direction only slightly. 
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Figure 12 shows the evolution of the parallel and perpendicular components of angular 
momentum J c /M 2 as a function of maximum density p max for models SFOO, SF45, and SF90. 
The evolution of these components is also similar to that for the moderate field models. It can 
be clearly seen that Jn and J± in model SF45 progress similarly to the respective components 
of SFOO and SF90, showing a strong selective effect of magnetic braking on J±. In all stages, 
J|l in model SFOO is 0.6 — 0.7 of that in model MF00 (shown for comparison), implying that 
the effect of magnetic braking on J\\ is greater under stronger magnetic fields, although the 
dependence is only minor. In fact, the dependence of the strength of the magnetic braking 
effect on field strength is more prominent for J±. For the strong field models, J±/M 2 also 
undulates during the runaway collapse, but the amplitude of fluctuation is considerably 
smaller than in the moderate field models. Furthermore, J±/M 2 in the strong field models 
decreases by a factor of ~ 100 prior to adiabatic core formation, whereas that in the moderate 
field models only decreases by a factor of ~ 10. 



5.3. Clouds Threaded by Weak Magnetic Fields 

The weak field models were initialized with a magnetic field of B = 7.42 /xG. Although 
the magnetic field of these models is very weak compared to the observed fields of dense 
molecular cloud cores (e.g., Crutcher 1999), these models show the dependence on the initial 
field strength. 

Figures 13 and 14 show the central cloud ejecting an outflow in the final stage for model 
WF45. The outflow has a bipolar bubble-shaped structure, distinguishable from that in 
model MF45 (Figs. 2c and 2d). In the outflow region, a toroidal magnetic field is prominent 
compared to the poloidal field, as shown in Figure 13a. This is another type of outflow driven 
by the magnetic pressure gradient of the toroidal field component (I-type flow in Tomisaka 
2002). The direction of outflow is no longer perpendicular to the disk of a 100 AU scale. 
The disk is considerably warped, and the normal direction of the disk varies with the radius 
of the disk. 

Figure 15 shows the warp of the disk quantitatively in terms of the normal vector of 
the disk p 3 with radius. A disk of 100 AU scale is inclined with respect to the z axis by 
40° toward the y direction, while a disk of a 10 AU scale is inclined by 20° in the opposite 
direction. The normal direction of the disk is roughly associated with the direction of the 
magnetic field, indicating that the magnetic field controls the warp of the disk. The direction 
of the normal vector of the disk is completely parallel to the magnetic field within 3 AU, 
corresponding to the radius of the adiabatic core. The direction of angular momentum also 
drifts in this plane, converging with that of the magnetic field within 3 AU. The slight change 
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in the direction of the outflow seen in Figure 14 seems to be driven by the change in the 
direction of the magnetic field outward from the adiabatic core. 

Figure 16 shows the loci of the magnetic field B c , the normal vector of the disk p 3 , 
and the angular momentum J c for a region of p > 0.1p max . The points indicating the three 
directions drift in this plane, with the loci of the magnetic field and the normal vector of the 
disk resembling those in Figure 15, suggesting that the configuration of the magnetic field 
and the density distribution reflect the history of cloud collapse. In the stage of adiabatic 
core formation, these directions are notably divergent. However, After outflow formation, 
the directions rapidly converge. 

This convergence of directions is also seen in Figure 17. For model WF45, J\\/M 2 
oscillates remarkably compared to J±/M 2 , attributable to the obtuse angle between J c and 
B c at 10~ 13 gcm -3 < p max < 10~ 12 gcm~ 3 . In the final stage, J« is considerably larger than 
J_i_, indicating that the angular momentum and magnetic field are aligned at this point. 
Figure 17 shows the evolution of J c in models WFOO and WF90. Compared with model 
MFOO, J\\/M 2 in model WFOO is much larger, indicating weaker magnetic braking in the 
weaker field model. 



6. Discussion 

6.1. Outflow Velocity 

The present simulations have shown that outflows form in the moderate and weak field 
models. Figure 18 shows the maximum outward radial velocity in the accretion phase for 
all moderate and weak field models. The evolution of the outflow velocity can be seen to 
depend on the field strength. In the moderate field models, small bumps occur at 150 yr < 
t — t c < 370 yr, representing bouncing in the central adiabatic core, where t c denotes the core 
formation epoch. The increase in the radial velocity in the interval t— t c > 370 yr corresponds 
to the period of maximum outflow velocity. The outflow velocity increases monotonically up 
until the final stage in models MFOO and MF45, but increases rapidly to 6 km s' 1 and then 
gradually decreases in model MF90. The maximum outflow velocities are roughly equal to 
the spin velocity of the outflows, implying that the outflow is accelerated by the magneto- 
centrifugal mechanism (c.f., Blandford & Payne 1982; Pudritz & Norman 1986). In contrast 
to the moderate field models, the outflow velocities in the weak field models exhibit little 
dependence on 9. The outflow velocities increase gradually immediately after adiabatic core 
formation, reaching maxima of ~ 2 kms -1 in the final stage, considerably slower than in the 
moderate field models. 
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Figure 19 shows the mean Alfven velocity in the outflow region. The mean Alfven 
velocity is defined by 

vA=1 vL.u,^ ir ' (2i) 

where V = f' ni dr denotes the volume of v r > 0.1c,. The Alfven velocities are almost 

JVr>0.1c s — 

identical in each field strength model. In the moderate models, the Alfven velocities increase 
rapidly at t — t c ~ 300 — 400 yr, synchronizing with the outflow formation, and decrease 
gradually to 3 — 4 km s^ 1 , which is approximately same as the outflow velocities in the final 
stages. The high-speed outflow shown in Figure 2d has velocity higher than the local Alfven 
velocity, while remaining region has sub-Alfvenic velocity. For the weak field models, the 
Alfven velocities are kept within 1 — 2 km s -1 after the outflow formation. This velocity is 
also comparable to the outflow velocity. 

The present simulations give maximum outflow velocities of ~ 6 km s _1 for the moderate 
field models, and ~ 2kms _1 for the weak field models. These velocities are significantly 
slower than the observed outflow velocities of ~ a few xlOOkms -1 (e.g., Richer et al. 2000; 
Bachiller & Tafalla 1999) considering the typical outflow speed of optical jets. Low-velocity 
outflow has also been reported in the numerical simulations of Tomisaka (2002, 1998) and 
more recently Allen et al. (2003). 

According to the steady solutions of a jet-disk system found by Kudoh & Shibata (1997), 
the outflow velocity increases in proportion to v 2 J^ . In the strong field models, Alfven 
speed increases up to the order of 100 km s _1 , i.e., 430c s = 82 km s _1 for model SF00, 
319c s = 61 km s^ 1 for model SF45, and 375c s = 71 km s _1 for model SF90 in the final stage, 
while the moderate field models give va = 3 — 4 km s -1 . Accordingly, the high-speed outflow 
of optical jets (~ 100 km s _1 ) are expected to occur in a cloud core with an initial magnetic 
field of ~ 40 fiG or more. 

Figure 20 shows the mechanical power of the outflow in terms of the radial momentum 
estimated by 



P = pv r dr , (22) 

J v r >5c s 

where v r denotes the radial velocity. Despite the slow outflow velocities in the weak field 
models (WF00 and WF45), the powers of the outflows are comparable to those of the mod- 
erate field models, implying that the volume of the outflow is also comparable to that of the 
moderate field models. The power is equal to ~ lO _4 M kms _1 in the final stage for model 
MF45, and appears to continue to increase. This increase in outflow power represents an 
increase in the volume or mass of the outflow region rather than an increase in the maximum 
velocity. The momentum flux is therefore Pj (t — t c ) ~ 1.6 x 10 -7 M km s _1 yr -1 in the final 
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stage for model MF45. However, this momentum flux is considerably smaller than the ob- 
served flux of > 10~ 6 M kms" 1 yr _1 for the bipolar molecular outflows around YSOs (see 
e.g., Bontemps et al. 1996). Tomisaka (2002) has shown that the momentum flux reaches 
10~ 6 — 10~ 5 M kms -1 yr -1 at t — t c ~ 10 4 yr, which represents much longer term evolution 
compared to the simulations in this study. 

We estimate the mass flux in the outflow (M out ) and the accretion rate (M in ) across 
the boundary of the 13th level grid, of which box-size is 71 AU, at the final stages. The 
outflow rates are M out = 6.56 x 10~ 7 M yr" 1 , 4.97 x 10~ 7 M yr" 1 , and 8.01 x 10~ 7 M yr" 1 
for models MF00, MF45, and MF90 while the accretion rates are identical to each other 
(Mi n = 4.1 x lO~ 5 M yr _1 ). A model with larger outflow rate shows larger outflow power 
P at the final stages. The ratio M out /M in ranges from 1.2% to 2.0%. The outflow rates in 
the weak field models are larger than those in the moderate field models such that M out = 
4.45 x 10~ 6 M yr" 1 , 3.90 x 10" 6 M yr" 1 , and 1.74 x 10" 6 M yr~ x for models WF00, WF45, 
and WF90. These are also in the same order of strength as the outflow power. The accretion 
rates are M[ n = (4.32 — 5.62) x lO _5 M yr _1 , which are approximately equal to that in 
the moderate field models. The ratio M out /M in ranges from 3.1% to 10%. The ratios of 
M out /M in presented here are considerably smaller than that of Tomisaka (2002) because of 
the early phase of outflow. 

6.2. Precession and Convergence of Magnetic Field, Angular Momentum, and 

Disk 

The cloud collapses with precession of the central magnetic field, angular momentum, 
and shape of the disk such that in the final stage, the directions of these descriptors converge 
to a single direction as the precession decreases in amplitude, as seen in Figures 5, 9, and 
16. The precession and convergence of the magnetic field and angular momentum is due 
to the interaction between them. Figure 21 shows a schematic diagram of this interaction, 
which consists of two mechanisms; rotation of the magnetic field and magnetic braking. 
The magnetic field threading the central cloud rotates due to the spin of the cloud, that 
is, the rotation axis of the magnetic field is aligned parallel to the angular momentum 
(Fig. 21a). The direction of the magnetic field B c drifts around a guiding center represented 
by the direction of the angular momentum «7 C , as seen in Figures 5, 9, and 16. The angular 
momentum is subsequently affected by magnetic braking. In the runaway collapse phase, 
J_i_ undulates as shown above (c.f., Mouschovias & Paleologou 1979), which also causes the 
direction of angular momentum to oscillate (Fig. 21b). In the outflow formation phase, the 
angular momentum rapidly inclines toward the local magnetic field, that is, the guiding 
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center moves due to B c , and the direction of J c follows, as shown in Figures 5, 9, and 
16. Consequently, the magnetic field and angular momentum precess until they come into 
alignment. 

The magnetic braking timescale can be estimated by a similar method to that employed 
in Mouschovias & Paleologou (1980) and Dorfi (1989). The Alfven wave transfers angular 
momentum from the central dense region to the outer less-dense region. Although the dense 
region of the collapsing magnetized cloud has a disk shape in the simulations, the central 
region of interest, i.e., the region of p > 0.1p max , is less oblate. Thus, the density distribution 
of the central collapsing cloud can be approximated by a spherical cloud for simplicity, as 
given by 

I Pc(r/\jY 2 forr^A, ' ^ 

where p c denotes the central density and Aj = c s (n / 'Gp c ) l l 2 denotes the Jeans length. The 
cloud has a central plateau and surrounding envelope, where the plateau maintains its scale 

— 1/2 

length to be roughly equal to the Jeans length, which decreases in proportion to p c ■ The 
z component of the moment of inertia of the plateau is estimated as 



I. = [ p(x 2 + y 2 )dr (24) 

Jr<Xj 

y^Pc^j ■ (25) 



Similarly, the moment of inertia of the envelope where the Alfven wave sweeps is estimated 

as 

J cxt (t) = f p(x 2 + y 2 )dr (26) 

J Xj<r<v a t 

3 



A 



^1 -1 



(27) 



where v a denotes the Alfven velocity. The braking time tb is estimated by equating the 
moments of inertia, 

h = Iext(t b ) , (28) 

leading to the braking time of 

The Alfven velocity v a remains roughly constant because the central magnetic field increases 
in proportion to pmax (see Fig. 3). However, the Jeans length Aj decreases in proportion to 
Pma( 2 in the runaway collapse phase, with the result that the braking time decreases during 
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the collapse in proportion to p ma £ . It should be noticed that the braking time has the 
same dependence p m a( 2 as that of the freefall timescale and the rotational period (without 
magnetic braking). 

The infalling envelope shows an infall velocity higher than the Alfven velocity in the 
runaway collapse phase, and one might image that the Alfven wave is trapped in the region 
of the central density plateau. The plateau shrinks, and its boundary moves inward at the 
speed higher than that of the outward propagating Alfven wave (v r + v a ). This indicates that 
the Alfven wave can escape from the central plateau and it transfers angular momentum to 
the infalling envelope. In the accretion phase, the Alfven velocity around the adiabatic core 
becomes extremely high (~ 100 km s -1 ), and the infalling gas has sub-Alfvenic velocity. 

Although the precession and convergence are attributed to nonlinear interaction be- 
tween magnetic braking and rotation, the magnetic braking during collapse can be roughly 
estimated by applying equation (29). Given the initial state of the present simulations, 
the braking times are t b = 1.64 x 10 5 yr, 3.28 x 10 6 yr, and 8.20 x 10 6 yr for the strong, 
moderate and weak field models, respectively. The precession is also related to the other 
timescales of the freefall time ts, and rotation timescale t TOt (= l/Q c )' tg = 2.10 x 10 5 yr and 
t rot = 1.28 x 10 6 yr in all models, where the rotation time is 6.6 times longer than the freefall 
time. In the strong field models, the braking time is considerably shorter than the freefall 
time, which allows magnetic braking to suppress the spin-up of the dense region during the 
collapse, as seen in Figure 11. In the moderate field models, the braking time is longer than 
the freefall time, resulting in considerable precession during the collapse as shown in Figures 
5 and 9. In the weak field models, the braking time is much longer than the freefall time, 
and the magnetic braking is too weak to redistribute the angular momentum significantly 
during the collapse. Moreover, the rotation timescale is comparable to the braking time, 
resulting in large precession, as shown in Figure 16. 

The precession of the disk axis is attributed to the precession of the magnetic field, which 
controls the collapse of the cloud by guiding the infalling gas. As the disk forms perpendicular 
to the direction of gas infall, the disk axis follows the direction of the magnetic field during 
the collapse. In this sense, the precession of the disk axis is passive. 

Although the cloud core is affected by magnetic braking in the simulations, the dense 
adiabatic core (the seed of a protostar) may be decoupled from the magnetic field as a result 
of efficient ambipolar diffusion (Nakano, Nishi, & Umebayashi 2002) in further stages. The 
directions of the disk, rotation axis, and outflow will be fixed at the convergence direction 
during the main accretion phase of the protostar. 
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6.3. Magnetic Braking and Torque 

The torque working on the dense region of p > 0.1p max during the collapse is estimated 
by performing f p>Q l r x (-)dr for equation (2), which yields 

£ = < 30 > 

and 

iVtot = N P + N m + N t + N 9 + N a . (31) 

The torques due to thermal pressure (Np), magnetic pressure (N m ), magnetic tension (N t ), 
and gravity (N g ), are then defined as 



N P = 

N m = 



I rx VPdr , (32) 

Jp>0.1p 

max 

/ r X V (^Pl dr ' ( 33 ) 

Jp>o.ip max V ° n J 

1 



N t = — (r x B)(B ■ dS) , (34) 

47T ,/p = o.lp max 

N g = - pr x V^rfr , (35) 

max 



r a = -/ p(rxv)(t>.dS). (36) 

./ P=0.1p max 



where f dS denotes the surface integral. The advection of angular momentum (JV ) is defined 
as 

N„ 

/p=0.1p max 

It should be pointed out that the volume of the above integral (region of p > 0.1p max ) 
changes temporally. In other words, the region of interest containing J c shrinks during 
the collapse. Accordingly, AT tot does not indicate an increasing rate of J c in its temporal 
evolution. However, it is true that AT tot represents the torque functioning on the region of 
interest, p > 0.1p max . 

Figure 22a shows the torques per mass |JV|/M as a function of p max for model MF45. 
The torque of magnetic tension N t and the advection of angular momentum iV a are much 
larger than the other components. The torque of magnetic tension N t causes magnetic 
braking, while the advection of angular momentum N a drives the spin-up of the central cloud. 
In the runaway collapse phase, N t /M remains almost constant, as indicated by equation (29) 
and the relationship of M oc pma( 2 - In the accretion phase, N t /M increases, promoting strong 
magnetic braking. Similarly, N a /M remains constant in the runaway collapse phase because 
vr oc R, j oc R, p oc R~ 2 , and M oc R. Therefore, M~ l J pjvpRdR is constant, where R 
and j denote the cylindrical radius and the specific angular momentum. 
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Figure 22b shows the parallel and perpendicular components of N t /M, N a /M, and 
N tot /M with respect to the angular momentum J c as a function of p max - The parallel and 
perpendicular components are defined as 

N\\ = \N-J C \/\J C \ , (37) 

and 

N± = l-ZV x J c |/| J c | . (38) 

The perpendicular component changes the direction of J c . In the runaway collapse phase, 
N a ^± is much smaller than iV a i|, indicating that gas accretion does not change the direction 
of rotation considerably. On the other hand, jV tj j_ remains roughly constant while N t » 
oscillates in the runaway collapse phase. In 10~ 16 gcm~ 3 < p max < 10~ 12 gcm -3 , A^ j_ is 
larger than N t> \\, with the result that the tension of the magnetic field changes the direction 
of J c considerably. 

In the early accretion phase, N tott ± is considerably larger than N tot u due to large iV tj j_. 
This torque drives J c to come into alignment with B c in the outflow formation epoch, 
resulting in the convergence of J c as shown in Figure 5. 



6.4. Application to Observations 

Observations of optical polarization reveal that the optical outflows tend to be aligned 
with the cloud-scale magnetic field. Cohen et al. (1984) showed that half of the well-defined 
outflows are aligned with the magnetic field within ±20°, and the remainder except for 
only one outflow are aligned within ±40° in the Cepheus A region. Similarly, Vrba et al. 
(1988) showed that the position angles of the outflows are inclined by 40° from the cloud- 
scale magnetic field in the Orion region. The results of the present simulations are highly 
consistent with these observations. 

The magnetic field lines inside the molecular cloud deviate significantly, as seen in 
Figures 7a and 15, and when the outflow extends out to > 100 AU, it is no longer parallel 
with the local magnetic field. The outflow velocity is comparable to the Alfven velocity, 
as seen in Figures 18 and 19. As the field strength is much weaker at the outer radius of 
the molecular cloud core, as seen in Figure 7b, the outflow can propagate straight in the 
direction in which the gas is launch at the root of the outflow. Accordingly, the outflow 
controls the direction of the magnetic field at the > 100 AU scale of the cloud core, as 
observed by Momose et al. (2001), although they argued that the magnetic field controls the 
outflow direction. 
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At the cloud-core scale, the configuration of the magnetic field is more complicated, as 
suggested by observations of submillimeter polarization (e.g., Wolf et al. 2003). In order 
to compare the present results with those observations, it is necessary to reproduce the 
polarization map using the simulation data. Such work is planned as part of future study. 

7. Summary 

The gravitational collapse of slowly rotating, magnetized molecular cloud cores was in- 
vestigated by high-resolution numerical simulations. Outflow formation in the collapsing 
cloud cores was also followed in the simulations. The cloud was found to collapse along 
magnetic field lines to form a disk perpendicular to the magnetic field. The central dense 
region spins up as a consequence of collapse. This spin-up is counteracted by magnetic 
braking, which was found to be more effective against the perpendicular component of an- 
gular momentum with respect to the local magnetic field than the parallel component. This 
leads to alignment of the magnetic field and rotation axis. With an initial magnetic field of 
~ 40 fiG, the magnetic field and rotation axis converge to alignment prior to the formation of 
an adiabatic core. With weaker initial magnetic fields, convergence also occurs rapidly due 
to the torque of magnetic tension associated with outflow formation. The central magnetic 
field is inclined from the cloud-scale magnetic field by no more than ~ 30° in the outflow 
formation stage given an initial magnetic field of ~ 20 fiG, and the inclination decreases 
with increasing initial field strength. 

The history of cloud collapse reproduces the distribution of magnetic field and density 
in the cloud core. The formation of a flat disk-like envelope at every scale of the cloud core 
occurs due to the perpendicular alignment of the disk with respect to the local magnetic 
field. In the envelope, the angular momentum is not required to be parallel to the local 
magnetic field, but becomes parallel at small scales of 10 AU. 

The outflow is always parallel to the magnetic field at the root of the outflow and 
perpendicular to the disk. In models with 9 < 80°, the central cloud converges to a near- 
axisymmetric model, and the outflow is attributed to the magneto-centrifugal wind. Even 
for models with 9 = 90°, the outflow is aligned parallel to the local magnetic field, and 
therefore perpendicular to the rotation axis. The acceleration is attributed to the pinch of 
the magnetic field. 

Consequently, disk-outflow systems are aligned with the cloud-scale magnetic field within 
~ 30° and ~ 5° for initial field strengths of ~ 20 fiG and ~ 40 /iG. This result is highly 
consistent with observations of the optical polarization of Cepheus A and Orion. 
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Fig. 1. — Evolution of moderate field model MF45. Density is shown in grayscale, and radial 
velocity is shown as arrows and red (positive) and blue (negative) contours at 0.2 km s _1 
intervals. The grid level is indicated by #Z. (a-c) Cross-sections taken along plane A at 
various time points, shown at different scales. (d,e) Cross-sections along (d) plane A and (e) 
plane B at the same time point. (/) Schematic showing planes A and B. 

Fig. 2. — Three-dimensional structures of density and magnetic field in the final stage of 
model MF45. Isosurfaces denote the iso-density surfaces of (a) p = p , (b) 10 2 po, and (c) 
10 5 po- Stream lines denote magnetic field lines, and the grid level is indicated by j^l. The 
box sizes are (a) 0.356 pc, (b) 4590 AU, and (c,d) 287 AU, where (d) shows the positive 
radial velocity (outflows) as v r = 10c s iso-velocity surfaces. 

Fig. 3. — Magnetic field strength and angular velocity in the dense region (p > 0.1p max ) 
as a function of maximum density p max f° r models MF00, MF45, and MF90. Upper solid 
lines denote magnetic field strength B c (c 2 po)~ 1 ^ 2 , the dotted line denotes the relationship 

1/2 

oc pmax, the dashed line denotes angular velocity of model NM, lower solid lines denote 
angular velocity f2 c (47rG r p )~ 1//2 , and diamonds denote the stage of adiabatic core formation 

(Pmax Per)- 

Fig. 4. — Angular momentum as a function of maximum density for models MF00, MF45, 
and MF90. Dashed line denotes J c /M 2 (c s /4irG) for model NM, solid lines denote angular 
momentum parallel to the local magnetic field J\\/M 2 (c s /4:7iG), dotted lines denote angular 
momentum perpendicular to the local magnetic field J±/M 2 (c s /4:TtG) and diamonds denote 
the stage of adiabatic core formation (p max = p cr )- 

Fig. 5. — Loci of directions of the central magnetic field B c (thin solid line), the central 
angular momentum J c (dotted line), the normal vector of the disk-like density structure 
p 3 (dashed line), and the magnetic field at the boundary B h (thick line) for model MF45. 
Diamonds denote the stage of adiabatic core formation (p max = p cr ), and triangles denote 
the final stage. (Inset) Detail of the convergence point. 

Fig. 6. — Loci of directions of the central magnetic field B c (thin solid line), the central 
angular momentum J c (dotted line), and the normal vector of the disk-like density structure 
p 3 (dashed line) with respect to the magnetic field at the boundary B^ for model MF45. 
Diamonds denote the stage of adiabatic core formation (p max = p cr ), and triangles denote 
the final stage. 
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Fig. 7. — Radial distribution of the mean magnetic field, mean angular momentum, and 
normal vector of the disk in the final stage of model MF45. (a) Directions of the magnetic 
field (solid line), angular momentum (dotted line), and normal vector of the disk (dashed 
line) for various radii. (6) Strength of mean magnetic field (solid line) and mean rotation 
velocity (dotted line) as a function of radius. 

Fig. 8. — Density, velocity and magnetic field distributions for model MF90. (a) Three- 
dimensional structure of the outflow region at t — t c = 464 yr, where the outflow speed 
reaches its maximum of 6.14 km s _1 . Streamlines denote the magnetic field lines, and the 
disk-like surface represents the iso-density surface of 10 5 po- Other surfaces denote the iso- 
velocity surfaces of v r = c s . Box size is 287 AU (/ = 12). (b) Cross-section of outflow region 
in the final stage using the same notation as Fig. 1. 

Fig. 9. — Loci of directions of the central magnetic field B c (thin solid line), the central 
angular momentum J c (dotted line), the normal vector of the disk-like density structure p 3 
(dashed line), and the magnetic field at the boundary B h (thick line) for models (a) MF70 
and (b) MF80. Diamonds denote the stage of adiabatic core formation (p max = p cr ), and 
triangles denote the final stage. 

Fig. 10. — Structure of strong field model SF45 in the final stage. Density is shown in 
grayscale, and radial velocity is shown as arrows and red (positive) and blue (negative) 
contours at 0.2 km s _1 intervals. The grid level is indicated by Cross-section is taken 
along plane B in Fig. If. 

Fig. 11. — Loci of directions of the central magnetic field B c (thin solid line), the central 
angular momentum J c (dotted line), the normal vector of the disk-like density structure 
p 3 (dashed line), and the magnetic field at the boundary By, (thick line) for model SF45. 
Diamonds denote the stage of adiabatic core formation (p max = Per), and triangles denote 
the final stage. 

Fig. 12. — Angular momentum as a function of maximum density for models SF00, SF45, 
and SF90. Dashed line denotes J c /M 2 (c s /4irG) for model NM, solid lines denote angular 
momentum parallel to the local magnetic field J\\/M 2 (c s /4irG), dotted lines denote angular 
momentum perpendicular to the local magnetic field J]_/M 2 (c s /4ttG), dot-dashed line de- 
notes J\\/M 2 for model MF00, and diamonds denote the stage of adiabatic core formation 

(Pmax Per)- 
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Fig. 13. — Three-dimensional structure of the outflow region in the final stage for model 
WF45. Box size is 287 AU (Z = 12). (a) Streamlines denote the magnetic field lines, and 
the disk-like surface represents the iso-density surface of 10 5 p - (b) Iso- velocity surface of 
v r c s . 

Fig. 14. — Structure of strong field model WF45 in the final stage. Density is shown in 
grayscale, and radial velocity is shown as arrows and red (positive) and blue (negative) 
contours at 0.2 km s _1 intervals. The grid level is indicated by #Z. Cross-section is taken 
along plane A in Fig. If. 

Fig. 15. — Radial distribution of the mean magnetic field (solid line), mean angular momen- 
tum (dotted line), and normal vector of the disk (dashed line) in the final stage of model 
WF45. 

Fig. 16. — Loci of directions of the central magnetic field B c (thin solid line), the central 
angular momentum J c (dotted line), the normal vector of the disk-like density structure 
p 3 (dashed line), and the magnetic field at the boundary B h (thick line) for model WF45. 
Diamonds denote the stage of adiabatic core formation (p max = Per), and triangles denote 
the final stage. 

Fig. 17. — Angular momentum as a function of maximum density for models WF00, WF45, 
and WF90. Dashed line denotes J c /M 2 (c s /47rG) for model NM, solid lines denote angular 
momentum parallel to the local magnetic field J\\/M 2 (c s /4irG), dotted lines denote angular 
momentum perpendicular to the local magnetic field Jj_/M 2 (c s /4nG), dot-dashed line de- 
notes J\\/M 2 for model MF00, and diamonds denote the stage of adiabatic core formation 

(Pmax Per)' 

Fig. 18. — Maximum radial velocities of outflows as a function of time from the adiabatic 
core formation epoch for models MF00, MF45, MF90, WF00, WF45, and WF90. 

Fig. 19. — Mean Alfven speeds in the outflow region (v r > 0.1c s ) as a function of time 
from the adiabatic core formation epoch for models MF00, MF45, MF90, WF00, WF45, and 
WF90. 

Fig. 20. — Radial momentum of the outflow region (y r > 5c s ) as a function of time from the 
adiabatic core formation epoch for models MF00, MF45, MF90, WF00, WF45, and WF90. 
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Fig. 21. — Schematic diagram of interaction between magnetic field and angular momentum, 
(a) The magnetic field rotates due to the spin of the cloud around the angular momentum 
vector. (6) The angular momentum inclines from the magnetic field and oscillates due to 
magnetic braking. 



Fig. 22. — Torque per mass N/M acting on the dense region of p > 0.1p max as a function 
of p max for model MF45. (a) Solid lines denote N t /M, N a /M, N m /M, N g /M, and N P /M, 
dashed line denotes N tot /M. (b) Thick, thin, and gray lines denote Nt/M, N a /M, and 
N tot /M. Solid and dashed lines denote parallel and perpendicular components with respect 
to angular momentum J c . 
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Table 1. Models and magnetic field strengths 



Models 


B (nG) 


Bo/Bcr 


Emag / 1 Egrav 


dv 


SFOO, SF45, SF90 


37.1 


0.50 


2.88 


0.658 


MFOO, MF45, MF70, MF80, MF90 


18.6 


0.25 


0.721 


2.63 


WFOO, WF45, WF90 


7.42 


0.10 


0.115 


16.5 


NM 


0.00 


0.00 


0.00 


oo 



